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ABSTRACT 

Cosmological simulations suggest that dark matter haloes are not spherical, but typically mod- 
erately to strongly triaxial systems. We investigate methods to convert spherical potential- 
density pairs into axisymmetric ones, in which the basic characteristics of the density profile 
(such as the slope at small and large radii) are retained. We achieve this goal by replacing the 
spherical radius r by an oblate radius m in the expression of the gravitational potential 3>(r). 

We extend and formalize the approach pioneered by Miyamoto & Nagai (1975) to be ap- 
plicable to arbitrary potential-density pairs. Unfortunately, an asymptotic study demonstrates 
that, at large radii, such models always show a R~ 3 disc superposed on a smooth roughly 
spherical density distribution. As a result, this recipe cannot be used to construct simple flat- 
tened potential-density pairs for dynamical systems such as dark matter haloes. Therefore 
we apply a modification of our original recipe that cures the problem of the discy behaviour. 
An asymptotic analysis now shows that the density distribution has the desired asymptotic 
behaviour at large radii (if the density falls less rapidly than r -4 ). We also show that the flat- 
tening procedure does not alter the shape of the density distribution at small radii: while the 
inner density contours are flattened, the slope of the density profile is unaltered. 

We apply this recipe to construct a set of flattened dark matter haloes based on the realistic 
spherical halo models by Dehnen & McLaughlin (2005). This example illustrates that the 
method works fine for modest flattening values, whereas stronger flattening values lead to 
peanut- shaped density distributions. 
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1 INTRODUCTION 



The intrinsic shape of dark matter haloes in the Cold Dark Matter paradigm has been a matter of debate for the past two decades (e.g. Frenk 
et al. 1988; Dubinski & Carlberg 1991; Warren et al. 1992; Cole & Lacey 1996). The most recent numerical simulations of cosmological 
clustering in a ACDM framework suggest that dark matter haloes are typically slowly rotating, moderately to strongly triaxial systems, with 
an average flattening that increases with increasing halo mass (e.g. Jing & Suto 2002; Hopkins, Bahcall & Bode 2005; Bailin & Steinmetz 
2005; Kasun & Evrard 2005; Allgood et al. 2006; Hayashi, Navarro & Springel 2007). In order to set up realistic models for dark matter 
haloes, we need to leave the well-trodden path of spherical symmetry and consider structures that reflect this variety in shapes. A crucial step 
is the construction of exact flattened or triaxial potential-density pairs with a realistic density structure (for example a cuspy behaviour at 
small radii). The most straightforward way to construct such axisymmetric or triaxial potential-density pairs is to consider a spherical density 
profile p(r) with the desired characteristics. If one replaces in this expression the spherical radius r by an oblate spheroidal radius 




(1) 



or a triaxial radius 




(2) 



one obtains an axisymmetric or triaxial density profile pf (r) = pf (m) whose isodensity surfaces are concentric oblate spheroids or ellipsoids 
with fixed axis ratios. This is a very straightforward recipe to construct simple oblate, prolate or triaxial density profiles, but unfortunately the 
formulae to calculate the corresponding potential <3>f (r) are cumbersome and can only be calculated numerically. In N-body or hydrodynam- 
ical simulations, the potential (or its gradient) has to be calculated very frequently and should be known with as efficiently and accurately as 
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possible to avoid numerical artifacts. As a result, it is worthwhile to look for methods to create exact potential-density pairs for axisymmetric 
or triaxial configurations in which the potential can be calculated in a more straightforward way. 

One general approach one can apply is to modify a spherical potential-density pair by adding a (finite) spherical harmonics expansion to 
it. While this approach can in principle be used to construct exact potential-density pairs for axisymmetric or triaxial systems (e.g. de Zeeuw 
& Carollo 1996), the number of terms in the expansion typically has to be substantial, such that we approach a similar situation as before. An 
alternative approach to construct axisymmetric or triaxial potential-density pairs is to apply a transformation r^mon the potential rather 
than on the density. If we have an analytical expression for <3>f(r) = <3>f(ra), we can use Poisson's equation to calculate the corresponding 
density distribution pf(r). The advantage of this approach is that the calculation of the density pf(r) from the potential ^f(ra) requires 
only differentiations, such that one always obtains potential-density pairs in which both potential and density are analytical functions. One 
disadvantage is that one cannot a priori set the shape of the isodensity contours; it is not even guaranteed that the density corresponding to 
a given potential is everywhere positive. A second disadvantage is that simple recipes such as equation (1) do not work for the potential, 
because the isopotential surfaces need to become spherical at large radii. Other, unavoidably more complicated, recipes must be considered. 

A successful example to construct an oblate analytical potential-density pair based on this strategy was the seminal work by Miyamoto 
& Nagai (1975). These authors create an exact potential-density pair by considering the spherical Plummer potential 

*^ = -v&< 

and turning it into a flattened potential 

MR, z) = — GM == . (4) 

\/R 2 + (a + Vz^+b 2 ) 

The resulting density distribution can be calculated fairly easily. This technique of constructing flattened or triaxial potential-density pairs has 
been applied to other potentials by e.g. Satoh (1980), Zhenglu (2000) and Vogt & Letelier (2007). While this recipe provides a straightforward 
means to create analytical flattened potential-density pairs, it generates densities that do not correspond to astrophysical systems such as dark 
matter haloes. Indeed, at large radii the Miyamoto-Nagai type models show a notorious double behaviour, with a thin disc superposed on 
a roughly spherical density profile. While this could provide a first approximation for disc galaxies (although the density in galaxy discs is 
exponential rather than R~ 3 ), this is obviously not the desired shape for a dark matter halo. 

The goal of the present paper is to construct axisymmetric or triaxial potential-density pairs that can be used to model realistic dark 
matter haloes, starting from an arbitrary spherical potential-density pair. In Section 2 we formulate a recipe based on the Miyamoto & Nagai 
(1975) approach and we demonstrate that such models always have the conspicuous R~ 3 disc-like behaviour at large radii, such that it 
cannot be used to represent dynamical systems as dark matter haloes. In Section 3 we adapt our recipe to cure this conspicuous behaviour 
and we demonstrate that this new recipe can produce analytical axisymmetric potential-density pairs which retain the original behaviour of 
the original density profile at both small and large radii. In Section 4 we apply this recipe to create a set of flattened potential-density suitable 
to represent dark matter haloes based on the spherical Dehnen & McLaughlin (2005) models. Section 5 sums up. 



2 MIYAMOTO-NAGAI TYPE MODELS 

In order to create a set of oblate spheroidal models, we will replace the spherical radius r in the spherical model by a spheroidal radius 
m = m(R, z). We need to find a prescript such that m is flattened at small radii and spherical at large radii. A general form that satisfies 
these requirements is 



rn(R, z) = |/i? 2 + g^y> (5) 

where Q(z) is a monotonically rising function of z which has the asymptotic behaviour that it approaches 1 at large radii and a constant 
value q at small radii. There are of course infinitely many different functions Q(z) with these characteristics. Inspired by Miyamoto & Nagai 
(1975) we concentrate on one such particular function 

Q(z) = = , (6) 



y/(a 2 -b 2 + Vz 2 + b 2 ) 2 



for which we find 



Substituting expression (6) into (5) we obtain 



m(R, z) = \Jr 2 + (a-b+^/z 2 + b 2 J 



(7) 



(8) 



When we move to polar coordinates in the equatorial plane using R — r sin and z — r cos 6, we find that at small radii 

V / cos 2 6' + ^ 2 sin 2 6> 

m(r, 0) ~ r, (9) 



Exact potential-density pairs for flattened dark haloes 3 



whereas at large radii m(r,Q) ~ r. The equipotential surfaces of the potential <£f (r, 6) — <3>f (m) will hence be flattened with flattening q at 
small radii and spherical at large radii, as required. 

Once we have constructed the potential <3>f (r, 0), we can compute the corresponding density pf(r, 0) that self-consistently generates this 
potential through Poisson's equation. 



1 



4 7rr 2 



<9r 



dm d$ f \ 1 5> 



<9r dm / sin 2 



. . 5m d$ f \ 



(10) 



After some manipulation we can write this equation as 

pt (r, 0) = Ti (r, 0) p{rn) + T 2 (r, 6) p(m) . (11) 

Here, p(m) represents the density of the original spherical model evaluated at r = m and p(m) is the average density of the original spherical 
model within a radius r, 

M(r) 



p(r) 



4ZL r 3 
3 1 



also evaluated at r = m. The coefficients Ti and T2 in equation (1 1) are independent of the potential of the system, and are given by 



and 



with 



T 2 (r,9) = 



3m 2 



Ti(r,0) 



2G 
F 2 



»+ 1 



TT 
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+ 
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cos 2 6> 
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These coefficients are never divergent and it is straightforward to check that for the trivial case a 
such that we recover the original spherical density distribution p(r), as required. 

The asymptotic expansion for the two coefficients T\ and T 2 at large radii reads 

b) |cos#| 



6, they reduce to T\ = 1 and T2 



Ti(r,t 



1 + 



r 2 (r, 0) 



r 



(12) 

(13) 

(14) 

(15) 
(16) 
(17) 
= 0, 

(18) 
(19) 



If the initial spherical model has a finite mass, the mean density p(m) behaves asymptotically as r -3 . If the density of the initial spherical 
model decreases more slowly than r -5 , the density of the corresponding oblate model will asymptotically reduce to the density of the 
spherical model. 

There is one caveat, however: the asymptotic expansion for the coefficients Ti and T2 derived above is not valid within the equatorial 
plane = f . Since z = in the entire equatorial plane, we have by construction F = q 2 and G = H = at every position R = r = m in 
the equatorial plane, such that we have the identities 



Ti (R, f ) = 1, (20) 

In the equatorial plane, the contribution of the second term in (1 1) will not decrease as R~ 5 , but only as R~ 3 , and it will therefore dominate 
the contribution of the first term. The mean density term p{m) is only dominant at large radii in the equatorial plane, such that the final density 
Pf (r, 6) will have a discy structure at large radii, with a thin R~ 3 disc superposed on a nearly spherical density profile. This behaviour is 
present in all Miyamoto-Nagai type models (Miyamoto & Nagai 1975; Satoh 1980; Zhenglu 2000) and unfortunately turns this recipe 
unsuitable for the construction of realistic flattened dark matter halo models. 



3 A NEW RECIPE TO CONSTRUCT FLATTENED MODELS 
3.1 Construction 

An inconvenience about the family of models presented in the previous section was the appearance of a discy structure at large radii. At large 
radii, the density profiles of the models behave similarly as their spherical progenitors, at all polar angles except in the equatorial plane z = 0. 
Mathematically, the reason for this discy structure is the discontinuous asymptotic behaviour of the function Q(z). Outside the equatorial 
plane, \z\ — > oo and Q — > 1 at large radii, such that the coefficient T2 disappears. In the equatorial plane, however, z is always equal to zero 
and T 2 does not converge to zero, such that the mean density term p(m) contributes to (and dominates) the density. 
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If we want to construct a set of models which do not show this discy structure at large radii, we can try to adapt the recipe from the 
previous section in a way that Q smoothly goes to one at large radii, both inside and outside the equatorial plane. This can be achieved by 
replacing the oblate radius m from equation (5) by 



m(R, z) 



R 2 + 



Using the functional form (6) for Q(r), we obtain 

m(R, z) 

or 

m(r, 0) - 



R 2 + 



(a-b + \/r 2 + b 2 y 



r 2 sin + cos 2 



(a-b+ \/r 2 + b 2 y 



(22) 



(23) 



(24) 



This new oblate radius m has a similar asymptotic behaviour at small and large radii as the original version (8). Replacing r by this oblate 
radius m in the expression <l>(r) for the potential of a spherical model, we obtain an oblate axisymmetric potential <3>f(r, 0). The density 
that self-consistently generates this potential can be found through Poisson's equation. Going through a similar analysis as in the previous 
Section, we find that the density can be written as 



^ (r, 0) = Ut (r, 0) p(m) + U 2 (r, 9) p{m) , 
where the coefficients Ui and U2 are now given by 



U 1 (r,9) = 



• 2 a . 1 — G 2 a , G 



1 , G \ *a 



and 



U 2 (r,0) = 



3 771 2 



1 + 1 + 



H 



G 



F 2 2F 2 F 3 



cos 2 ^+^ 



H G 2 H G 2 \ 

1 , 

2 F 2F 4F 2 J 



with now obviously 

G = rF'(r), 
H = r 2 F"(r). 

Again, it is elementary to check that for Ui = 1 and U2 = for the trivial case a = 6. 



Q>) 



(25) 
(26) 

(27) 

(28) 
(29) 
(30) 



3.2 The asymptotic behaviour at large radii 

After some algebra, one obtains that the asymptotic expansion of the coefficients Ui and U2 at large radii becomes 

Ui(r,0)~l, (31) 
2(a-|)co^ 

These expansions are valid for all directions; in particular, there is no different behaviour of these components in the equatorial plane. 
Comparing the asymptotic behaviour of these coefficients to those of the previous Section, we see that the coefficients T\ and U± both have 
an asymptotic behaviour with simply 1 as the leading term, such that the first term in the total density (11) or (25) behaves identically as 
the density of the original spherical model. The coefficients T2 and U2 behave in a different way, however: while T2 decreases as r~ 2 at 
large radii, U2 decreases only as r _1 , with an additional cos 20 azimuthal dependence. Only for those models where the density of the initial 
spherical model decreases more slowly than r -4 , the density of the corresponding flattened model will asymptotically reduce to the density 
of the spherical model. For models where the density decreases more rapidly at large radii, the density will be dominated by the second 
term. As this term becomes negative at large radii near the rotation axis, such models are unacceptable. This method to flatten spherical 
potential-density pairs is therefore limited to those models where the density decreases more slowly than r -4 . This is typically the case for 
dark matter haloes. 



3.3 The asymptotic behaviour at small radii 

The goal of this paper was to consider a method to, starting from a spherical potential-density pair, construct a flattened analog that preserves 
the asymptotic behaviour of the density profile at both small and large radii. An asymptotic analysis of the functions U\ and U2 learns that 
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both terms converge to a finite value at small radii, 



1 cos + q sin 6 

Ui(r,0) g — - — - . 2 , (33) 

<T cos 2 + <? 2 sin 2 



-r-r / J- 17 U Dill 1/ KsVJQ 1/ ^„ 

f/ 2 (r, ) ~ -5-f- — . (34) 



1-q 2 q 2 sin 2 6> - 2 cos 2 
Sq 2 cos 2 + q 2 sin 2 
Assume that our density profile p(r) behaves as a power law at small radii, 

p(r)~p r- 70 , (35) 
with < 70 < 3. It is easy to check that the mean density p(r) has a similar slope, 

(36) 

O—70 

Combining the expressions (9), (25), (33), (34), (35) and (36), we find the asymptotic behaviour of pf (r, 0) at small radii, 

n (r f)) p0 (! - 7o^ 2 + 2g 2 ) Q 2 sin 2 + (1 - 70 + 1q 2 ) cos 2 

W ' } ~ 3 - 70 4 2 -™ (cos 2 + q 2 sin 2 0)i+7o/2 ' V> 

This result demonstrates that the asymptotic behaviour of the original spherical model is preserved after the flattening: the slope of the density 
profile remains, only the zero point changes and becomes dependent of the polar angle 0. In particular, along the major and minor axes we 
find the expansions 

^-Cl^ *"" t38) 

Using these expressions, we easily find the expression for the flattening q p of the isodensity surfaces at small radii, 

\1 - 7o<T + 2q 2 J 
For modest flattening (q < 1) we find a nearly linear relation, 

Apart from the well-known fact that the isodensity surfaces are generally more flattened than the isopotential surfaces, this equation illustrates 
that the flattening of the inner isodensity surfaces increases with increasing 70 . 



4 APPLICATION: FLATTENED DEHNEN-MCLAUGHLIN DARK HALO MODELS 

As an application of the formalism we have described in the previous section, we construct in this section a set of flattened dark matter halo 
models. Our starting point is the general three-parameter family of generalized NFW models or Zhao models 

i+ yj • (42) 

The density of this general model is characterized by a power law with (negative) slopes 70 and 700 at respectively small and large radii. 
The third parameter, 5, is a measure of the width of the transition region between the two power-law zones. This three-parameter model, 
the general dynamical properties of which are described by Saha (1993) and Zhao (1996), includes a large variety of well-known simple 
potential density pairs used extensively throughout the literature to describe dynamical systems ranging from galactic nuclei to dark matter 
haloes (e.g. Plummer 1911; Jaffe 1983; Hernquist 1990; Dehnen 1993; Tremaine et al. 1994; Navarro, Frenk, & White 1997). 

To demonstrate our recipe to construct oblate dark matter haloes, we focus on a particular, non-trivial subset of this general three- 
parameter set of models presented by Dehnen & McLaughlin (2005). These authors used the Jeans equations to derived the density profile of a 
self-gravitating spherically symmetric dynamical system that satisfy two characteristics of dark haloes observed in cosmological simulations. 
The first of these characteristics is the mysterious power-law behaviour of the so-called pseudo phase-space density 1 Q — p/a 3 with a(r) 
the velocity dispersion (Taylor & Navarro 2001; Ascasibar et al. 2004; Rasia, Tormen, & Moscardini 2004). The second ingredient in their 
model is a linear relation between the logarithmic density slope and the anisotropy, observed in different haloes with radically different 
formation histories (Hansen & Moore 2006; Hansen & Stadel 2006). Constrained by these two assumptions, together with the assumption of 



pir) = po 



1 In fact, Dehnen & McLaughlin (2005) assumed a power-law behaviour of the quantity p/ 1 a% with a r the radial velocity dispersion and e a free parameter. 
We assume e = 3, the most natural choice. 
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Figure 1. Isopotential (top row) and isodensity (bottom row) plots in the meriodional plane of the Dehnen & McLaughlin (2005) halo models flattened 
according to the Miyamoto-Nagai type flattening from Section 2. Different models are shown according to different values of the potential flattening parameter 
q, ranging from 1 to 0.7. Each plot has dimensions of 10a x 10a. At large radii, the isodensity surfaces are discy with an R~ 3 dependence in the equatorial 
plane. 



spherical symmetry, Dehnen & McLaughlin (2005) derived a dark matter halo density of the form (42), with 

7 + 10A) 



7o 



31 - 20A) 
9 

4-2$, 



po = 



9 ' 
1 - 5/3 M 
18tt r% 



(43) 
(44) 
(45) 
(46) 



with (3o the value of the velocity anisotropy parameter at small radii and M the total mass. The gravitational potential that generates this 
mass density is 



$(r) 



GM 



B 



11 - 10 A) 



(47) 



4-2/3 r i+(r/ rQ )v \4-2f3o 1 4 - 2(3 

with B x (p, q) the incomplete Beta function. 

We constructed flattened versions of the spherical Dehnen & McLaughlin (2005) haloes with j3o — and ro = a using the recipes of 
Sections 2 and 3. In Figure 1 we plot the isopotential and isodensity surfaces in the meridional plane for the Miyamoto-Nagai type models for 
different values of the parameter q = b/a. The isopotential surfaces (top row) are nearly spheroidal, and the flattening decreases smoothly 
from the centre (where the axis ratio is equal to q) to spherical at large radii. The bottom row shows the isodensity surface plots of the 
corresponding models. The isodensity surfaces are nearly spheroidal in the central regions with a flattening that is much stronger than the 
flattening of the potential. For increasing radii, the isodensity surfaces become increasingly lemon-shaped, and at large radii, they become 
very discy with an R~ 3 dependence in the equatorial plane, in agreement with the results from Section 2. Obviously, these models can not 
represent realistic dark haloes. 

In Figure 2 we plot in a similar way the isopotential and isodensity surfaces for flattened haloes according to the recipe of Section 3. 
Compared to Figure 2, the shape of the isopotential surfaces is hardly different, with just a small alteration (a smoother behaviour) near the 
equatorial plane. This minor change does affect the large-scale behaviour of the isodensity surfaces significantly: the discy structure at large 
radii disappears and the isodensity surfaces remain roughly spheroidal with a flattening that smoothly disappears for increasing radii. 

For models with a modest potential flattening (0.9 < q < 1), all isodensity surfaces are nearly spheroidal, with a flattening parameter 
q p ranging from roughly 0.8 at the centre to 1 at large radii. For stronger flattening parameters, the shape of the isodensity surfaces at small 
radii becomes increasingly more peanut- shaped, as progressively more mass has be located near the equatorial plane (and progressively less 
mass near the symmetry axis) in order to stretch the isopotential surfaces into an oblate shape. This net transfer of mass from the symmetry 
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Figure 2. Similar to Figure 1, but now the models are flattened according to the adapted recipe from Section 3. The shape of the isopotential surfaces is slightly 
smoother near the equatorial plane, which causes the discy structure at large radii in the isodensity surfaces to disappear. Instead, the isodensity surfaces remain 
roughly spheroidal with a flattening that smoothly disappears for increasing radii. 



axis to the equatorial plane cannot continue for ever. At a certain stage, a negative density will arise around on either side of the equatorial 
plane, causing physically unacceptable models. 



5 DISCUSSION AND CONCLUSION 

In this paper, we have investigated methods to convert spherical potential-density pairs into axisymmetric ones, in which the basic charac- 
teristics of the density profile (such as the slope at small and large radii) are retained. We attempted this by replacing the spherical radius 
r by an oblate radius m in the expression of the gravitational potential <3>(r). The advantage of this approach is that the calculation of the 
corresponding density via Poisson's equation requires only differentiations, such that one always obtains fully analytical potential-density 
pairs. Disadvantages are that one cannot a priori set the shape of the isodensity surfaces and that the standard recipe for the oblate radius m 
cannot be used since the isopotential surfaces need to become spherical at large radii. 

In Section 2 we have considered a recipe inspired by the flattening of the Plummer potential-density pair by Miyamoto & Nagai (1975). 
We have extended and formalized this mechanism to be applicable to arbitrary potential-density pairs. The recipe is sufficiently simple: 
the density can formally be written as a simple linear combination of the original spherical density p(m) and the mean density p(m) - 
both evaluated at the oblate radius m- and the coefficients are simple, non-divergent functions independent of the potential-density pair. 
Unfortunately, an asymptotic study demonstrates that, at large radii, such models always show a R~ 3 disc superposed on a smooth roughly 
spherical density distribution. As a result, this recipe cannot be used to construct simple flattened potential-density pairs for dynamical 
systems such as early-type galaxies or dark matter haloes. 

In Section 3 we have applied a modification of our original recipe that cures the problem of the discy behaviour. The new models can be 
constructed in a similar manner as the Miyamoto-Nagai type models from Section 2; only the coefficients in the linear combination of p(m) 
and p(m) are different. An asymptotic analysis now shows that the density distribution has the desired asymptotic behaviour at large radii (if 
the density falls less rapidly than r -4 ). We also show that the flattening procedure does not alter the shape of the density distribution at small 
radii: while the inner density surfaces are flattened, the slope of the density profile is unaltered. We have applied this recipe to construct a set 
of flattened dark matter haloes based on the realistic spherical halo models by Dehnen & McLaughlin (2005). This example illustrates that 
the method works fine for modest flattening values (0.9 < q ^ 1 with q the flattening of the isopotential surfaces at small radii), whereas 
stronger flattening values lead to peanut- shaped density distributions. 

For the sake of simplicity, we have focused this work on the construction of flattened potential-density pairs, but there is no reason why 
the recipe applied here should be limited to an oblate geometry. The procedure is readily applicable to prolate or triaxial geometries, which 
seems necessary to represent the general population of dark matter haloes. 
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